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Abstract. It is proven that a class of Gegenbauer tau approximations to a 4th order differential 
eigenvalue problem of hydrodynamic type provide real, negative, and distinct eigenvalues, as is the 
case for the exact solutions. This class of Gegenbauer tau methods includes Chebyshev and Legendre 
Galerkin and 'inviscid' Galerkin but does not include Chebyshev and Legendre tau. Rigorous and 
numerical results show that the results are sharp: positive or complex eigenvalues arise outside of 
this class. The widely used modified tau approach is proved to be equivalent to the Galerkin method. 
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1. Introduction. The Chebyshev tau method used by Orszag [T7] to obtain 
exponentially accurate solutions of the Orr-Sommerfeld equation yields two eigenval- 
ues with large positive real parts. Such eigenvalues also occur for Stokes modes in a 
channel given by the fourth order differential equation 



with the boundary conditions u — Du = at x = ±1, where A is the eigenvalue, 
u = u(x) is the eigenf unction, D = d/dx and a is a real wavenumber. The Stokes 
eigenvalues A arc real and negative as can be checked by multiplying (II. lj) by u* , the 
complex conjugate of u{x), and integrating by parts twice using the no-slip boundary 
conditions. In fact the Stokes spectrum has been known analytically since Rayleigh [JJ 
§26.1]. Yet, the Chebyshev tau method applied to (jl.ip yields 2 eigenvalues with large 
positive real parts, for any order of approximation and for any numerical accuracy. 
Such eigenvalues are obviously spurious for Gottlieb and Orszag Chap. 13] 

introduced the eigenvalue problem 



as an even simpler ID model of incompressible fluid flow. This is the a — > limit of 
the eigenvalue problem 1)) and of the Orr-Sommerfeld equation [17] . For any fixed 
a, problem (| 1 . 2[) is also the asymptotic equation for large A solutions of the Stokes 
and Orr-Sommerfeld equations. The eigensolutions of (|1.2| are known analytically. 
They consist of even modes u{x) = 1 — cos(n7rx)/ cos(n7r) with A = — n 2 ii 2 and 
odd modes u(x) = x — sin(g n x)/sin(g„), with A = —q„, where q n = tan<?„, so that 
mr < q n < (2n + l)7r/2, V integer n > 0. The key properties of these solutions 
are that the eigenvalues are real, negative, and distinct, and the even and odd mode 
eigenvalues interlace. These properties also hold for Stokes modes, the solutions of 
(jl.ljl . but not for Orr-Sommerfeld modes. 

The Chebyshev tau method provides spectrally accurate approximations to the 
lower magnitude eigenvalues but it also yields two large positive eigenvalues for prob- 
lems (jl.ip and (|1.2p [TTl Table 13.1]. Those positive eigenvalues are clearly spurious 
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(1.1) 



{D 2 ~a 2 fu = X(D 2 -a 2 )u 



(1.2) 



D 4 u = \D 2 u in - 1 < x < 1, 
u = Du = at x = ±1, 
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Fig. 1.1. Eigenvalues of a 3D steady state solution of the Navier-Stokes equations for plane 
Couette flow \21^j computed with Chebyshev tau (left) and Chebyshev Galerkin (right) for identical 
resolutions (8773 modes after symmetry reductions. The solutions themselves are indistinguishable) . 
Note the difference in horizontal scales. Chebyshev tau produces 274 eigenvalues with positive real 
parts, 273 of which are spurious. Chebyshev Galerkin returns only one positive eigenvalue, the 
physical one, equal to 0.03681 at Reynolds number 1000 1211 Fig. 4]- 



since it is known that (jl.ip and (|1.2|) should only have negative eigenvalues. The 
Chebyshev tau method yields 2 spurious eigenvalues for no-slip (a.k.a. 'clamped') 
boundary conditions u = Du = at x = ±1, but none for the free-slip boundary 
conditions u — D 2 u — at x = ±1. The latter problem reduces to the 2nd order 
problem D 2 v = Xv with v(±l) = for which a class of Jacobi and Gegenbauer tau 
methods has been proven to yield real, negative and distinct eigenvalues H] . For 
mixed boundary conditions, e.g. u(±l) = Du{— 1) = D 2 u(l) = there is one spurious 
eigenvalue (this is a numerical observation). 

For a ID problem such as (jl.lj) or (|1.2|) . the spurious eigenvalues are easy to 
recognize and they appear as minor nuisances. Boyd [5J §7.6] even questions the value 
of distinguishing between 'spurious' and numerically inaccurate eigenvalues. However, 
in many applications, large negative eigenvalues are inconsequential, while 'spurious' 
positive eigenvalues are very significant, and in higher dimensions, spurious eigenvalues 
are not as easy to pick out and set aside. In a recent application, 3D unstable traveling 
wave solutions of the Navier-Stokes equations were calculated with both free-slip and 
no-slip boundary conditions, and anything in between, by Newton's method [191 120] . 
In that application, the Chebyshev tau method provides hundreds of spurious unstable 
eigenvalues, depending on resolution and the exact type of boundary conditions, not 
all of which have very large magnitudes (fig. left). A simple change in the test 
functions from Chebyshev polynomials T n (x) to (1 — x 2 )T n (x) or (1 — x 2 ) 2 T n (x) 
eliminates all those spurious eigenvalues (fig. 11.11 right). This, and more, is proven 
below for the test problem (|1.2[) in the broader context of Gegenbauer tau methods 
which include Chebyshev and Legendre tau, and Chebyshev and Legendre Galerkin 
methods. 

Another practical consequence of the spurious eigenvalues is that the Chebyshev 
tau method is unconditionally unstable when applied to the time-dependent version 
of (jl.lj) or (|1.2p . with d/dt in place of A. Such time-dependent problems appear 
as building blocks in Navier-Stokes simulations of channel-type flows. Gottlieb and 
Orszag [111 p. 145] proposed a modified tau method for the time-dependent problems 
and proved that the modified method was stable for even solutions. The modified 
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tau method (sect. [5]) is a key idea behind several successful time integration schemes 
for the Navier-Stokes equations [5J §7.3], [T3]. The modified tau method amounts 
to using 2 more expansion polynomials for the 4th order differential operator on the 
left hand sides of (| 1 . 2[) and than for the 2nd order operator on the right hand 
sides. That modified tau method was adapted to eigenvalue problems by Gardner, 
Trogdon and Douglass [8 j and McFadden, Murray and Boisvert [16j . McFadden et al. 
showed the equivalence between the modified Chebyshev tau method and a Chebyshev 
Galerkin method by direct calculation. Zcbib [H] had given numerical evidence that 
the Galerkin method removed spurious eigenvalues. The modified tau method idea 
was adapted to the collocation formulation by Huang and Sloan [15] . 

A heuristic 'explanation' for spurious eigenvalues is that there is a 'mismatch' 
between the number of boundary conditions applied to the 4th order operator on the 
left hand side of (|1 .2|) and those applied to the 2nd order operator on the right hand 
side. That interpretation fits with the modified tau method which uses two more 
polynomials for the 4th order operator than for the 2nd order operator. However it is 
incorrect since, while the tau method for Chebyshev polynomials of the 1st kind T n (x) 
gives spurious eigenvalues, for instance, the tau method for Chebyshev polynomials 
of the 2nd kind U n (x) does not. 

All of these various methods are best seen in the context of the Gegenbauer 
class with residuals weighted by W^'(x) = (1 — x 2 ) 1 ^ 1 / 2 (sect. [2]), where 7 = 
corresponds to Chebyshev and 7 = 1/2 to Legendre polynomials. Dawkins et al.]§\ 
proved existence of spurious positive eigenvalues for (jl.2[) when 7 < 1/2. The proof 
is straightforward. For (|1.2p . the polynomial equation for fi = 1/A can be derived 
explicitly (sect. [4|. All coefficients of that polynomial are real and positive, except 
the constant term which is negative when 7 < 1/2. Hence there is one real positive 
[i, and a 'spurious' positive eigenvalue, when 7 < 1/2 (details are given in sect. 
15. 2p . For 7 = 1/2, the Legendre tau case, the constant term is zero, hence there is 
one /i = eigenvalue, or a A = l//i = 00 eigenvalue. Perturbation analysis shows 
that the A = 00 eigenvalues become very large positive eigenvalues for 7 < 1/2 
and very large negative eigenvalues for 7 > 1/2. We provide a quicker derivation 
of those results in section [3] Dawkins et al.'s results do not prove that there are 
no spurious eigenvalues for 7 > 1/2 since there could be complex eigenvalues with 
positive real parts. In section [5[ we prove that the Gegenbauer tau method applied to 
(II. 2|) provides eigenvalues that are real, negative, and distinct when 1/2 < 7 < 7/2. 
This provides a complete characterization of the Gegenbauer tau spectrum for problem 
(11. 2| . Numerical calculations confirm that the range 1/2 < 7 < 7/2 is sharp. Spurious 
positive eigenvalues exist for 7 < 1/2 [F] and complex eigenvalues arise for 7 > 7/2, 
for sufficiently high polynomial order. In section [6[ we prove that the modified tau 
method is mathematically equivalent to the Galerkin approach. 

Obviously, 7 = 1/2 is a critical value for the weig ht function (1 - x 2 )^- 1 ! 2 . The 
boundaries x — ±1 have infinite weight for 7 < 1/2, and zero weight for 7 > 1/2, 
but we do not know a valid heuristic explanation for 'spurious' eigenvalues beyond 
that observation, if one exists. Section [3J provides further insights into the nature 
of the spurious eigenvalues and gives some support for the view that 'spurious' and 
numerically inaccurate eigenvalues are related. In figure fTTTl for fixed resolution, the 
273 spurious eigenvalues for 7 = (Chebyshev tau), escape to +00 as 7 / 1/2 
(Legendre tau). They come back from —00 as 7 increases beyond 1/2. Thus, there is 
indeed a connection between large positive eigenvalues and large negative eigenvalues, 
but whether we have 'spurious' positive eigenvalues or inconsequential very negative 
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eigenvalues is sharply controlled by 7, irrespective of the order of approximation n 
(see also eqn. (|3.6[1 ). 

2. Tau and Galerkin methods. Definition 2.1. A Gegenbauer tau method 
approximates the solution u(x) of a differential equation in — 1 < x < 1 by a polyno- 
mial of degree n, u n (x), that satisfies the m boundary conditions exactly. The remain- 
ing n + 1 — m polynomial coefficients are determined by imposing that the residual be 
orthogonal to all polynomials of degree n — m (or less ) with respect to the Gegenbauer 
weight W^\x) = (1 - 2; 2 ) 7 - 1 / 2 , with 7 > -1/2. 

For problem (|1.2[) , the residual 

(2.1) Rn-2(x) = XD 2 U n (x) - D 4 U n ( X ) 

is a polynomial of degree n — 2 in x. The polynomial approximation u n (x) is deter- 
mined from the 4 boundary conditions u n (±l) = Du n (±l) — and the requirement 
that R n _ 2 {x) is orthogonal to all polynomials q n —4,( x ) °f degree n — 4 or less with 
respect to the weight function W^\x) = (1 — a; 2 ) 7-1 / 2 > in the interval (— 1, 1) 

(2.2) J R n - 2 (x)q n - A (x)W^(x)dx = 0, Vq n - 4 (x). 

This provides n — 3 equations which together with the 4 boundary conditions yield 
n+ 1 equations for the n+1 undetermined coefficients in the polynomial approximation 
u n (x). For the Gegenbauer weight function (x) = (1 — x 2 ) 7 ^ 1 / 2 , the residual can 
be written explicitly as 

(2.3) R n _ 2 (x) = r Q XG { :l 2 (x) + nXG^x) 

for some x- independent coefficients tq and t\ , where G^ (x) is the Gegenbauer poly- 
nomial of degree n. This follows from orthogonality of the Gegenbauer polynomials 
in — 1 < x < 1 with respect to the weight (1 — x 2 ) 7-1 / 2 which implies orthogonality 
of the Gegenbauer polynomial of degree k to any polynomial of degree k — 1 or less 
with respect to that weight function. 

Gegenbauer (a.k.a. ultraspherical) polynomials are a special subclass of the Ja- 
cobi polynomials pQ. The latter are the most general class of polynomial solutions 
of a Sturm-Liouville eigenproblem that is singular at ±1, as required for faster than 
algebraic convergence [3] . Gegenbauer polynomials are the most general class of poly- 
nomials with the odd-even symmetry G^?\x) = (— l) n Gn?\— x). This is a one param- 
eter family of polynomials, with the parameter 7 > —1/2. Chebyshev polynomials 
correspond to 7 = and Legendre polynomials to 7 = 1/2. We use a (slightly) 
non-standard normalization of Gegenbauer polynomials since the standard normal- 
ization [1] is singular in the Chebyshev case. Some key properties of Gegenbauer 
polynomials used in this paper are given in Appendix El Note that if A = then, 
from (|2.ip . the residual R n -2(x) must be a polynomial of degree n — 4 implying that 
To = T\ = in (|2.3p and D A u n {x) = for all x in (—1, 1). The boundary conditions 
w n (±l) = Du n (±l) = then imply that u n {x) = for all x in [—1,1], the trivial 
solution. Hence we can assume that A ^ in the Gegenbauer tau method applied to 
& 

Following common usage (Tl] [3l [22l [16] , we have 

Definition 2.2. A Gegenbauer Galerkin method approximates the solution u(x) 
of a differential equation in —1 < x < I by a polynomial of degree n, u n (x), that 
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satisfies the m boundary conditions exactly. The remaining n+l — m polynomial co- 
efficients are determined by imposing that the residual be orthogonal to all polynomials 
of degree n (or less ) that satisfy the homogeneous boundary conditions, with respect 
to the Gegenbauer weight W^^x) = (1 - a; 2 ) 7 " 1 / 2 , with 7 > -1/2. 

Strictly speaking, this a Petrov-Galerkin method since the test functions are not 
identical to the trial functions because of the Gegenbauer weight (1 — a; 2 ) 7-1 / 2 [3]. 

For problem (| 1 .2(1 . u n {x) is determined from the boundary conditions u n (±l) = 
Du n (±l) = and orthogonality, with respect to weight W^ J \x) = (1 — a: 2 ) 7-1 / 2 , 
of the residual (|2.1[) to all polynomials of degree n that vanish together with their 
derivative at x = ±1. Such polynomials can be written as (1 — x 2 ) 2 q n -n{x) where 
q n -i(x) is an arbitrary polynomial of degree n— 4, and the weighted residual equations 
read 

(2.4) J R n ^ 2 (x)(l-x 2 ) 2 q n ^(x)W^(x)dx = 0, Vq n ^(x). 

The Gegenbauer Galerkin method is therefore equivalent to the tau method for the 
weight W^ +2 \x) = (1 — x 2 ) 2 W^ 7 )(a;) and its residual has the explicit form 

(2.5) R n . 2 {x)=T XG^\x)+r 1 \G^{x). 

So a Chebyshev (or Legendre) Galerkin method for clamped boundary conditions, 
w n (±l) = Du n (±l) = 0, is in fact a tau method for Chebyshev (or Legendre) polyno- 
mials of the 3rd kind (proportional to the 2nd derivative of Chebyshev (or Legendre) 
polynomials (|A.5[) V Since we consider a range of the Gegenbauer parameter 7, the 
Gegenbauer tau method also includes some Gegenbauer Galerkin methods. 

This suggests an intermediate method where the test functions are polynomials 
that vanish at x = ±1 (inviscid boundary conditions only). 

Definition 2.3. The Gegenbauer 'inviscid Galerkin' method determines u n (x) 
from the 4 boundary conditions u„(±l) = Du n (±l) = and orthogonality of the 
residual to all polynomials of degree n — 2 that vanish at x = ±1, with respect to the 
weight function W^\x) = (1 — a; 2 ) 7-1 / 2 . 

Such test polynomials can be written in the form (1 — x 2 )q n -4 L (x) where qn-^x) 
is an arbitrary polynomial of degree n — 4, so the weighted residual equations read 

(2.6) J R n - 2 (x)(l-x 2 )q n _ 4 (x)W^(x)dx = 0, Vq n -4(x). 

The Gegenbauer inviscid Galerkin method is therefore equivalent to a Gegenbauer 
tau method with weight W^-' 1+1 \x) and its residual for (II. 2|) is 

(2.7) i?„-2 (x) = to A (x) + ri A (x). 

Thus, a Chebyshev (or Legendre) inviscid Galerkin method is a tau method for Cheby- 
shev (or Legendre) polynomials of the 2nd kind. Since we consider a range of the 
Gegenbauer parameter 7, the Gegenbauer tau method also includes some Gegenbauer 
Inviscid Galerkin methods. 

For completeness, we list the collocation approach, where u n (x) is determined 
from the boundary conditions u„(±l) = Du n (±l) = and enforcing i?„_2(^j) = at 
the n — 3 interior Gauss-Lobatto points Xj such that DG n -%(xj) = 0, j = 1, n — 3, 
[3 §2.2]. The residual has the form [TU1 eqn. (4.5)] 



(2.8) 



R n - 2 {x) = (A + Bx) DG n - 2 (x), 



G 
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for some A and B independent of x. That residual can be written in several equivalent 
forms by using the properties of Gegenbauer polynomials (appendix IX)) . We do not 
have rigorous results for the collocation method. 

3. Legendre and near-Legendre tau cases. Here we provide a quicker and 
more complete derivation of earlier results [6] about spurious eigenvalues for the Legen- 
dre and near-Legendre tau case. This section provides a useful technical introduction 
to the problem but is not necessary to derive the main results of this paper. Dawkins 
et a/.[B], focusing only on even modes, use the monomial basis x 2k to derive an explicit 
form for the generalized eigenvalue problem Aa = XBa for the Legendre tau method. 
In the monomial basis, the matrix A is upper triangular and nonsingular, and the 
matrix B is upper Hessenberg but its first row is identically zero, hence there exists 
one infinite eigenvalue. A perturbation analysis is used to show that the infinite eigen- 
value of the Legendre tau method becomes a large positive eigenvalue for Gegenbauer 
tau methods with 7 < 1/2 and a large negative eigenvalue for 7 > 1/2. 

In the Legendre tau method, the polynomial approximation u n (x) of degree n 
to problem (|1.2|) satisfies the 4 boundary conditions u n — Du n — at x = ±1. 
Thus u n (x) = (1 — x 2 ) 2 p n -i{x) and the polynomial p n _ 4 (a;) is determined from the 
weighted residual equations (|2.2p with 7=1/2 and W^^ 2 '(x) = 1, 



(3.1) 



(/J.D u n - D 2 u n ) q n ~i{x)dx — 0, \/q n -i{x). 



The mathematical problem is fully specified, except for an arbitrary multiplicative 
constant for u n {x). Choosing various polynomial bases for p n -i(x) and q n -i{x) will 
lead to distinct matrix problems but those problems are all similar to each other and 
provide exactly the same eigenvalues, in exact arithmetic. 

We use the bases G\ f \x) for p n -^(x) and G^ 2 \x) for q n -i(x), with k,l — 
0, . . . ,n — 4, where G^ > \x) is the Gegenbauer polynomial of degree n for index 7 
(see appendix IA1 and recall that G^ 2 \x) = Pk{x) are Legendre polynomials). Thus 

we write u n (x) — J2"=o a i (1 ~ x 2 ) 2 G t f^ 2 \x), for some n — 3 coefficients ai to be 
determined. The tau equations (|3.1[) provide the matrix eigenproblem /iAa = Ba, or 

A 1 Ya=o A ( k > l ) a l = J2\ l =o B ( k > l ) a l witn 



(3.2) 
(3.3) 

for k,l = 0,. 

(3.4) 

(3.5) 



A(kJ) 



D 



B(k,l) = I D 



-1 
1 



{l-x 2 ) 2 Gf ,2 \x) 



(l-x 2 ) 2 Gf /2 \x) 



G^ /2) {x) dx, 



G { t /2) {x) dx, 



4. Using (|B.1[) . these expressions simplify to 



A(k, l)=Cij \p 2 G { ^ 2) {x)\ G[ 1/2) (x) dx, 
B(k,l)=d [ G^ 2 \x) G[ 1/2 \x) dx, 



where d = i(i + 1)(Z + 2){l + + 4). 

P n {x) are orthogonal with respect 
yields that B(k,l) cx ^^+2, where 5k,i+2 is the 



Since the Legendre polynomials G„^ 2 \x) 



to the unit weight, equation 
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Kronecker delta, so that 5(0, 1) = 5(1, 1) = for all I and 5 has non zero elements 
only on the second subdiagonal. For A(k,l), use (| A. 10|) to express D 2 G^^ ( x ) as a 
linear combination of G, 1 (x) , G^_ 2 f (x), etc. Orthogonality of the Legendre poly- 
nomials Gn (x) then implies that A(k, I) is upper triangular with non-zero diagonal 
elements. Hence A is non-singular while the nullspace of 5 is two-dimensional. The 
eigenvalue problem [iAa = Ba therefore has two (i — eigenvalues. Since the only non 
zero elements of B consists of the sub-diagonal B(l + 2, 1), the two right eigenvectors 
corresponding to \i — are a = [0, . . . , 0, 1] T and [0, . . . , 0, 1, 0] T . In other words, 

(3.6) u n {x) = {l-x 2 ) 2 G^{x) and u n (x) = (1 - x 2 ) 2 G^lf (x) 

satisfy the boundary conditions and the tau equations (|3.1[) with fi = = 1/A, for 
all n > 5. One mode is even, the other one is odd. Likewise, the left eigenvectors 
b T = [1, . . . , 0] and [0, 1, . . . , 0], satisfy fib T A = b T B with /i = 0. These results are 
for the Legendre case and /i — corresponds toA = l//x=oo. 

Now consider the Gegenbauer tau equations for 7— 1/2 = e with |e| <C 1, the near- 
Legendre case. The equations arc (|3.ip but with the extra weig ht factor W {l) {x) = 
(1 — x 2 ) e inside the integral. We can figure out what happens to the /x = eigenvalues 
of the e = Legendre case by perturbation. The matrices A and B and the left 
and right eigenvectors, denoted a and b respectively, as well as the eigenvalue /i now 
depend on e. Let 

(3.7) A = Ao + eAi+Oie 2 ), B = B + eB l + 0(e 2 ), 

(3.8) a = a + e ai +0{e 2 ), b = b + eb x + 0(e 2 ), 

(3.9) a* = Mo + +0(e 2 ), 



where A$ and Bq are the matrices obtained above in (13. 4ft and (13. 5|) for e = while 
do and by are the corresponding right and left eigenvectors so that /loAoao — Bqclq 
and noby Ao = b^Bo. Substituting these e-expansions in the eigenvalue equation 
fiAa — Ba and canceling out the zeroth order term, we obtain 

(3.10) HiA a + noAxao + ^oA ai = Bxa + B ax + 0(e). 

Multiplying by by cancels out the ^b^ Aqcl\ = b$ Bqcli terms so we obtain 

b^Bicio - (j, bffAia 



(3.11) Ml 



b^A a 



This expression is general but simplifies further since we are interested in the pertur- 
bation of the zero eigenvalues fj,o — 0- This expression for /ii is quite simple since bo, 
ao and Aq are the zero" 1 order objects. All we need to compute when /iq = is the 
first order correction B\ to the matrix B. But since ao and 60 have only one non-zero 
component as given at the end of the previous paragraph, we only need to calculate 
two components of the B matrix. For n even, all that is needed are the first order 
corrections to B(0,n — 4) for the even mode and to B(l,n — 5) for the odd mode. 
For n odd, we need B(0, n — 5) for the even mode and 5(1, n — 4) for the odd mode, 
however since even and odd modes decouple in this problem, it suffices to compute 
both even and odd modes in only one case of n even or odd. The matrix elements in 
the e ^ cases are still given by (I3.4p and (|3.5p but with the extra (1 — x 2 ) e weight 
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factor inside the integrals. Since G^' 2 '(x) = x and Gn^ 2 \x) — P n (x), the Legendre 
polynomial of degree n, we obtain 

(3.12) B(Q, n - 4) =C„_ 4 J P n -2{x) (1 - x 2 f dx = eB x (0, n - 4) + 0(e 2 ), 

(3.13) B(l,n-5)=C n - 5 J x P n - 3 {x) (1 - x 2 f dx = eB^l, n - 5) + 0(e 2 ), 

with C; as defined in (|3.5p . The integrals are readily evaluated and details are provided 
in appendix [Bj Using (|B.2|) . (|B.5|) and (IB.8|) we obtain for the even mode (for n even) 
that 

(3.14) Ml -^- 4 )- " 4 



Ao(0,n-4) (n-2) 2 (n- l) 2 ' 

This matches the formula in Dawkins et al. [6l page 456] since their 2N = n — 4 and 
2v — 1 = 2e. Likewise, using (|B.7j) . and (|B.9j) for the odd mode (with n even) yields 

B x (l,n-5) -4 

(3.15) jtii — 



A (l,n-5) (n-4) 2 (n- l) 2 ' 

Again, if n is odd then /ii for the even mode is given by (|3.14p but with n — 1 in 
lieu of n. Likewise for n odd, the odd mode is given by (|3.15p with n + 1 in lieu of 
n. Finally, since A = l//i, the A = oo eigenvalues in the Legendre tau case, become 
A = l/(e/ii + 0(e 2 )) - l/(e/ii) in the near-Legendre cases. From (|3.14[) and (|3 . 1 5|) . 
these eigenvalues will be 0(n 4 /e). Furthermore they will be positive when e < 
(i.e. spurious when 7 < 1/2) but negative when e > 0. 

4. Characteristic Polynomials. For the model problem (|1.2|) , we can bypass 
the matrix eigenproblem of section [3] to directly derive the characteristic polynomial 
for the eigenvalues /j, = 1/A. To do so, invert equation (|2.1|) to express the polynomial 
approximation D 2 u n (x) in terms of the residual R n ^2{x) 

00 

(4.1) D 2 u n (x) = nJ2» kD2kR ™-2( x ) 

k=0 

where /i = 1/A. The inversion (|4.1|) follows from application of the geometric (Neu- 
mann) series for (1 — /iD 2 )^ 1 = J^'kLo ^ k D 2k which terminates since R n -2{x) is a 
polynomial. Thus, u n (x) can be computed in terms of the unknown tau coefficients 
by double integration of (|4.1[) and application of the boundary conditions. We can as- 
sume that A ^ because A = with u n (±l) = Du n (±l) = necessarily corresponds 
to the trivial solution u n (x) = 0, Vx in [—1,1], as noted in the previous section. 

The Gegenbauer polynomials are even in x for n even and odd for n odd (|A.8|) . 
The symmetry of the differential equation (|1.2[) and of the Gegenbauer polynomials 
allows decoupling of the discrete problem into even and odd solutions. This parity 
reduction leads to simpler residuals and simpler forms for the corresponding charac- 
teristic polynomials. The residual in the parity-separated Gegenbauer case contains 
only one term 

(4.2) R n ^(x)=T \G%l 2 (x), 
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instead of (|2.3[) . where Gn(x) is the Gegenbauer polynomial of degree n with n 
even for even solutions and odd for odd solutions. Substituting (|4.2|) in (14. ip and 
renormalizing u n (x) by tq gives 



(4.3) D^u n (x) = J2^D 2k G^l 2 (x). 



,k r>2fc/c(T) 
k=Q 

For 7 > 1/2, the identity in the form 2~fG%l 2 (x) = DG^~^\x) can be used to 

write (14.31) in the form 



oo 

(4.4) D*u n (x) = ±-Y j ^D^G^\x), 

~ fc=0 

which integrates to 

oo 

(4.5) Du^ar) = -^^^Glr/'W +C, 

~ fc=0 

where C is an arbitrary constant. 

4.1. Even Solutions. For even solutions u n (x) = u n (—x), n is even and Du n (x) 
is odd so C = in (|4.5[) . The boundary condition £>u„(l) = gives the characteristic 
equation for \i (for n even and 7 > 1/2) 



(4.6) ^ M fe ^ 2fe Gir i 1) (l) = 

fc=0 



4.2. Odd Solutions. For odd solutions, u n (x) = —u n (—x), n is odd and the 
boundary condition Du n (l) = requires that 

00 

(4.7) C = ~ y £^D M G^ 1 1 \l). 

' fe=0 

Substituting this C value in (|4. 5|) and integrating gives 

00 oc 

(4.8) 2 7M „(x) = E^Z^Gfci 1 ^) - x£fLP($-V(l) 

k=0 k=0 

where we must define 

(4.9) = f Gfc 1 ^-) *• = Gir y~ Gt of iX) 

Jo 2(n + 7-2) 

since u n (x) and n are odd, where we have used (|A.10|) to evaluate the integral and 
the symmetry (|A.8j) so that G„(0) = G„_2(0) = for n odd. The boundary condition 
u n (l) = yields the characteristic polynomial equation (for n odd and 7 > 1/2) 

oc 00 

(4.10) ^^G^U) -^ k D 2k G^\l) = 0. 

fe=0 fe=0 
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For 7 > 3/2, we can use identity in the form 2( 7 -l)G, ( 1 I 1 1) (x) = L>G.[, 7 2) (x), 
to write the characteristic equation (|4.10p as 

oc oo 

(4.11) J2^ D2kG n' 2) ( 1 )-T / ^ E>2k+lG n~ 2) ^)=0- 
k=0 k=0 

For 1/2 < 7 < 3/2, this cannot be used since 7 — 2 < —1/2, but using (|4. 9|) for the 
D^ 1 term in the first sum, the characteristic equation (|4. 1 0[) can be written 

(4.12) ^ kD2k+lG %-?\ i )- \<:\ \ -E^^w = °- 

5. Zeros of Characteristic Polynomials. Here we prove that the zeros of the 
characteristic polynomial equations (14. 6p and (|4.10[) are real, negative, and distinct 
for 1/2 < 7 < 7/2. Some background material is needed. 

5.1. Stable Polynomials and the Hermite Biehler Theorem. A polyno- 
mial p(z) is stable if and only if all its zeros have negative real parts. Stable poly- 
nomials can arise as characteristic polynomials of a numerical method applied to a 
differential equation as in [5] for Du — Xu with u(l) = and in other dynamical 
systems applications. The characterization of stable polynomials that is most useful 
here is given by [H], [HI p. 197], 

Theorem 5.1. The Hermite-Biehler Theorem. The real polynomial p(z) = 
f2(z 2 ) + z9(z 2 ) is stable if and only i/f2(/i) and 0(/i) form a positive pair. 

Definition 5.2. Two real polynomials f2(/it) and 0(/i) of degree n and n — 1 (or 
n ) respectively, form a positive pair if: 

(a) the roots /ii,/i 2 , ■ ■ ■ , M„ of fi(^) and fi[, n' 2 , . . . , //^ (or fi[, /j,' 2 , . . . , fj! n ) of 
0(/i) are all distinct, real and negative. 

(b) the roots interlace as follows: fXi < fi± < /j,2 < • • • < /U _ x < < ( or 

fj>i < fii < " • < fj, < fi n < ) 

(c) the highest coefficients of £l(fj,) and 0(/i) are of like sign. 

We will use the following theorem about positive pairs [HI pl98], 
Theorem 5.3. Any nontrivial real linear combination of two polynomials of 
degree n (or n and n — 1) with interlacing roots has real roots. (Since such a linear 
combination changes sign n—1 times along the real axis, it has n— 1 real roots. Since 
it is a real polynomial of degree n, the remaining root is real also.) 

5.2. Eigenvalues for even modes. In [5] and [4] we study the Gegenbauer tau 
method for D 2 u — Xu with u(±l) = 0, which leads to the characteristic polynomials 
J2kLo f- k D 2k G^\ l) . The derivation of that result is entirely similar to that in sections 
[2] and [4] The strategy to prove that the Gegenbauer tau method for that 2nd order 
problem has real, negative and distinct roots is to show stability of the polynomial 

n 

(5.1) p(z) = J2 zkDkG n ) (l) 

k=0 

for —1/2 < 7 < 3/2 then to use the Hermite Biehler Theorem to deduce that 
Theorem 5.4. For — 1/2 < 7 < 3/2, the polynomials 

oc 00 

(5.2) n£V) = X>^ 2fc <4 7) (l) ™d 6«( M ) = £ ^ k D 2k+1 G^\l) 

k=0 k=0 



GEGENBAUER TAU METHODS WITH AND WITHOUT SPURIOUS EIGENVALUES 11 



form a positive pair. From HA.5\) this is equivalent to stating that the polynomials 

oo oo 

(5.3) nW( M )=£;/D a *GW(l) and = £ ^D 2k G^\l) 

also form a positive pair. Combining the 7 and 7+1 ranges in i f 5. 3\) yields that 
has real, negative, and distinct roots for —1/2 < 7 < 5/2. Stability of (|5.ip is proven 
in [5J Theorem 1] for the broader class of Jacobi polynomials. The basic ideas of the 
proof are along the lines of Gottlieb [|TJ and Gottlieb and Lustman's work [TIT] on 
stability of the Chebyshev collocation method for the 1st and 2nd order operator. 

For the 4th order problem (JT2J), D A u = XD 2 u with u(±l) = Du(±l) = 0, the 
Gegenbauer tau method gives the characteristic polynomial (|4.6[) for even solutions. 
This is the polynomial Q^Si\fJ-) of (|5.3p that appears for the 2nd order problem [JEJE] 
and is known to have real, negative, and distinct eigenvalues for —1/2 < 7— 1 < 5/2, 
that is for 1/2 < 7 < 7/2. Hence, it follows directly from (|4.6p and theorem 15.41 that 
the Gegenbauer tau approximation for even solutions of problem (|1.2I) . D A u — XD 2 u 
with u(±l) = Du(±l) — 0, has real, negative, and distinct eigenvalues for 1/2 < 7 < 
7/2. 

This result is sharp. For 7 > 7/2 and sufficiently large n, our numerical com- 
putations show that the polynomial has a pair of complex eigenvalues. For 7 < 1/2 
the polynomial (|4.6[) has a real positive eigenvalue as first proven in [5|. The proof 
goes as follows. For 7 < 1/2, we cannot use (|4.5[) since 7 — 1 is below the range 
of definition of Gegenbauer polynomials (appendix 1X1) . Instead, integrate (|4. 3|) and 
use identity (lA~T0l) to obtain G ( ^_ 2 dx = (G^ x (l) - G^ 3 ( l))/(2 (7 + n - 2)) since 
G„_i(0) = G„-3(0) = for n even (|A.8|) . Using formula (|A.11[) , one shows that 
increases with n if 7 > 1/2 but decreases with n if 7 < 1/2. Hence, the con- 
stant term is negative for 7 < 1/2 while all the other coefficients of the characteristic 
polynomial can be shown to be positive using (|A.5[) and (|A.11[) . Therefore, there is 
one real positive eigenvalue as proved in [6j. For 7 = 1/2, the Legendre case, the 
constant term vanishes and there is a /1 = (A = 00) eigenvalue as established in 
section [31 

Remark 1. Exact even solutions of D A u — XD 2 u with u(±l) = Du(±l) = 
obey D 3 u = XDu + C with C = 0, since D 3 u and Du are odd. Thus even solution 
eigenvalues of 11.2(1 are equal to the eigenvalues for odd solutions of the 2nd order 
problem D 2 w = Xw with w = at x — ±1, with w — Du. The Gegenbauer tau 
version of this property is that eigenvalues for even Gegenbauer tau solutions of il.2] ) 
of order n (even) and index 7 are equal to the eigenvalues for odd Gegenbauer tau 
solutions of D 2 w = Xw, w(±l) = of order n—1 (odd) and index 7— 1. This follows 
directly from |^.g[ ) and \J$. 

5.3. Eigenvalues for odd modes. The reduction of the 4th order problem (|1.2p 
to the 2nd order problem does not hold for odd modes which have the characteristic 
equation (|4.10[) . In fact, all previous theoretical work focused only on the even modes 
[51 [TT] . The same general strategy as for the even case led us to prove the stability of 
a shifted version of polynomial (|5.1[) . 

Theorem 5.5. Let Gn (x) denote the non standard Gegenbauer polynomial of 
degree n as defined in avvendix PA then the polynomial 

(5.4) *W W = b!=lW G »+ + £ z * D ^\l) 

Z[n + 1 > fe=0 
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is stable for -1/2 < 7 < 1/2. 

The proof is elementary but technical, it is given in Appendix[C] We also need the 
following simple lemma. This lemma will help us determine the sign of the coefficients 
of the characteristic polynomials. 

Lemma 5.6. With Gn\x) as defined in avvendix the expression 

(5.5) D k+1 G^(l)-D k G^(l) > 

for k = 0, . . . , n — 1 and 7 > —1/2. 
Proof. From (jJT3|) 

1 J " l J (n-fc)!r( 7 +l)r(27 + 2fe)' 

where r(z) is the standard gamma function. For k > and given that 7 > —1/2, the 
sign of the above expression is positive since all individual terms are positive. In the 
k = case, the sign of the expression is determined by the terms pP_ since all other 
terms are positive. If —1/2 < 7 < both numerator and denominator are negative 
and thus their ratio is positive. If 7 > the two terms are positive and thus again 
their ratio is positive. For 7 = 0, a simple limiting argument shows positiveness again, 



in fact from (TAll) . g£ 0) (1) = T n (l)/n = 1/ 



n. 



Taking the ratio D k+1 G^ ] (1) / D k Gi 7) (1) and making some simplifications gives 

f5 7) D k+1 G^\l) = (2 7 + n + fc)(n-fc) 

' j i>Gl 7) (l) ~ (2 7 + 2/c + l) 

Since k < n - 1 then 27 + 2k + 1 < 27 + (n - 1) + k + 1 = 27 + n + k. Thus 

(5.8) > 1 fc = 0...n-l 

£ fc Gi 7) (i) 

and since both derivatives are positive, the lemma follows. □ 
We now have all the tools to prove 

Theorem 5.7. The Gegenbauer tau approximation to problem hl.2\) has real, 
negative, and distinct eigenvalues for 1/2 < 7 < 7/2. This 7 range is sharp, spurious 
positive eigenvalues exist for 7 < 1/2 and complex eigenvalues arise for 7/2 < 7. 

Proof. This has already been proven in section 15.21 for even solutions. For odd 
solutions, we need to consider two separate cases. 

Case 1. 3/2<7<7/2. The characteristic polynomial (|4.11[) 

00 00 

(5.9) ^^^Gi^flj-^^^G^W ^^(^-e^ 2 )^) 

k=0 k=0 

is a linear combination of the polynomials f^ 7 ~ 2 '(/i) and Qn~ 2 \p), which form a 
positive pair for 3/2 < 7 < 7/2, by theorem 15.41 Therefore, by theorem 15.31 this 
characteristic polynomial has real roots. Then by lemma [B~6l we deduce that all its 
coefficients are of the same sign, hence all its roots must be negative. 
Case 2. 1/2 < 7 < 3/2. The polynomial 

(5.10) P { n-x\ z ) = A (^ 2 ) + z ®(z 2 ) 
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with p"n \z) as in theorem l5.51 is stable for the desired range of parameters by theorem 
15.51 so the Hermite Biehler theorem (theorem 15. ip implies that the polynomials 



( 5 - n ) am = Gil t ( l!+ 7 G y )(1) +E^^ 2fc gi 7 :i 1) (i) 



fe=0 



2(n + 7-2) 



and 

oo 

(5.12) =J2^D 2k+1 G^\l) = Q%?(p) 

fc=0 

form a positive pair, with f2^(/x) and Qn / \fi) as defined in theorem l5.41 Thus M$(/x) 
and A(/x) have interlacing roots and any real linear combination of the two must have 
real roots (theorem l5.3|) . Now the characteristic polynomial (14. 12)) is in fact the linear 
combination fJ,$>((A) — A(/i) so it has real roots. Its constant term is equal to 

(5 ' 13) 2(n + 7 -2) (1) 

which is negative if -1/2 < 7- 1 < 1/2, that is 1/2 < 7 < 3/2, from (|A~TTj) . All other 
coefficients of the characteristic polynomial — A(/z) are negative by lemma 

Since all coefficients have the same sign and all roots are real, all the roots must be 
negative. □ 



6. Galerkin and Modified tau methods. Our main theorem 15.71 can be ex- 
pressed in terms of the inviscid Galerkin and Galerkin methods since these methods 
are equivalent to Gegenbauer tau methods with index 7 + 1 and 7 + 2, respectively, 
as shown in section [5J 

Corollary 6.1. The Gegenbauer Inviscid Galerkin approximation to problem 
kl.Sfy has real, negative, and distinct eigenvalues for —1/2 < 7 < 5/2. 

COROLLARY 6.2. The Gegenbauer Galerkin approximation to problem tl.ty) has 
real negative eigenvalues for — 1/2 < 7 < 3/2. 

COROLLARY 6.3. Since Chebyshev corresponds to 7 = and Legendre to 7 = 1/2, 
the Chebyshev and Legendre tau approximations to problem have spurious eigen- 

values, but the Chebyshev and Legendre inviscid Galerkin (7=1 and 3/2, respectively) 
and Galerkin (7 = 2 and 5/2, respectively) approximations provide real, negative, and 
distinct eigenvalues. 

Finally, we prove that the modified tau method introduced by Gottlieb and Orszag 
pTj and developed by various authors [8], [16] is equivalent to the Galerkin method. 
McFadden, Murray and Boisvert [16j have already shown the equivalence between the 
modified Chebyshev tau and the Chebyshev Galerkin methods. Our simpler proof 
generalizes their results to the Gegenbauer class of approximations. 

The idea for the modified tau method, widely used for time-marching, starts with 
the substitution v(x) — D 2 u(x). Problem (|1.2p reads 



(6.1) 



D 2 u — v, D 2 v = Aw, with u = Du = at x = ±1. 
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If we approximate u(x) by a polynomial of degree n, then v = D 2 u suggests that v 
should be a polynomial of degree n— 2, however the modified tau method approximates 
both u(x) and v(x) by polynomials of degree n, 



(6.2) 



, n {x) =J2u k G^(x), D 2 u n (x) = £4 2) 4 7) (^ 

n n—2 

o n (x) = j2ti k G%\x), D 2 v n (x) = J2€ )G k ) (x), 



k=0 k=0 



where the superscripts indicate the Gegenbauer coefficients of the corresponding 
derivatives. These can be expressed in terms of the Gegenbauer coefficients of the 
original function using (|A.10[) twice, as in the Chebyshev tau method [3], |17) . Hence 
there are 2n + 2 coefficients to be determined, uq, . . . , u n , vq, . . . , v n . In the modified 
tau method, these are determined by the 4 boundary conditions and the 2n — 2 tau 
equations obtained from orthogonalizing the residuals of both equations D 2 u = v 
and D 2 v — Xv to the first n — 1 Gegenbauer polynomials G^\x), . . . , G^}_ 2 {x) with 
respect to the Gegenbauer weight (1 — x 2 ) 7 ~ 1//2 . In terms, of the expansions (|6.2[) . 
these weighted residual equations have the simple form 

u (2) =i) k < k <n-2, 
(6.3) - - 

v ( k > =\v k , 0<k<n-2. 

McFadden et al.\T<o\ showed that the modified Chebyshev tau is equivalent to the 
Chebyshev Galerkin method for this particular problem. We provide a simpler proof 
for the general setting of Gegenbauer polynomials. 

Theorem 6.4. The modified Gegenbauer tau method proposed in fTJjj is equivalent 
to the Gegenbauer Galerkin method for problem 

Proof. Let the polynomial approximations and their derivatives as in (|6.2[) . Then 
the tau equations (|6.3[) are equivalent to the residual equations 

Vn(x) - D 2 u n {x) = Vn^G^x) + v n G™(x), 
(X-D^Vn^^Xv^G^ix) +Xv n G^(x). 

Combining the two yields 
(6.5) 

(A - D 2 )D 2 u n (x) = D 2 G^\{x) + v n D 2 G^\x), u n (±l) = Du n {±l) = 
which using (|A.5|) is equivalent to 

(6.6) (A - D 2 )D 2 u n (x) = t XG [ ^\x) + nXG^ (x), u n (±l) = Du n {±\) = 0. 

with i> n = 4A(7 + 1)(7 + 2)tq and = 4A(7 + 1)(7 + 2)ti. This is exactly the 

Gegenbauer Galerkin method, as given in equation (|2.5p . □ 

Since the modified Gegenbauer tau is equivalent to the Galerkin method, the 
results in section [5] imply that the modified tau method for problem 11.21 has real and 
negative eigenvalues for —1/2 < 7 < 3/2. This includes Chebyshev for 7 = and 
Legendre for 7 = 1/2. 
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Appendix A. Gegenbauer (Ultraspherical) Polynomials. 

The Gegenbauer (a.k.a. Ultraspherical) polynomials Cn\x), 7 > —1/2, of degree 
n are the Jacobi polynomials with a — (3 = 7 — 1/2, up to normalization [TJ 22.5.20]. 
They are symmetric (even for n even and odd for n odd) orthogonal polynomials 
with weig ht function W^\x) = (1 - z 2 )^ 1 / 2 . Since the standard normalization [TJ 
22.3.4], is singular for the Chebyshev case 7 = 0, we use a non-standard normalization 
that includes the Chebyshev case but preserves the simplicity of the Gegenbauer 
recurrences. Set 



(A.l) 



C<?\x) 
2 7 : 



n > 1. 



These Gegenbauer polynomials satisfy the orthogonality relationship 



(A.2) 



where [TJ 22.2.3], 



(A.3) 



h 1 



tt2- 2 T(27+ 1) 



and 



hi 



0, m 7^ n, 
h%, m = n, 



7r2 -i-2 7r( - n + 2 7 ) 



r 2 ( 7 + l) ' n (n + 7 )n!r 2 ( 7 +l) 

The Sturm-Liouville form of the Gegenbauer equation for Gn^ (x) is 



n > 1. 



(A.4) 



D 



(1 - x 2 )^ +1 / 2 DG^\x)} = -n{n + 2 7 ) (1 - x 2 y~ 1/2 G^\x) 



They satisfy the derivative recurrence formula 

d 



(A.5) 



dx 



(?W =2( T + l)Gh+ 1) , 



(for C„ this is formula [3J A.57]), and their three-term recurrence takes the simple 
form 

(A.6) (n+l)G^ 1 = 2(n + 7 )^ ) -(ri-l + 27)4-i ) " > 2 - 

with 



(A.7) 



G™(x) = l, G["(x)=x, G^ = ( 7 + 1).t 2 - 



(7) 



1 



These recurrences can be used to verify the odd-even symmetry of Gegenbauer poly- 
nomials [TJ 22.4.2] 



(A.; 



G^(x) = (-l) n G^(-x). 



Differentiating the recurrence (|A.6|) with respect to x and subtracting from the cor- 
responding recurrence for 7 + 1 using (|A.5[) . yields [TJ 22.7.23] 



(A.9) 



(n + 7 )<#> = (7 + l) [^-^ 



n > 3. 
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Combined with (|A.5[) . this leads to the important derivative recurrence between 
Gegenbauer polynomials of same index 7 



(A.10) 



G$»\x) = DG<{>\x), 2(1 + 7 )G^(.t) = DG™(x), 



2(n + 7 )GfW 



d r 

dx 



G 



(7) 

n+l 



G 



(7) 

n-l 



Evaluating the Gegenbauer polynomial at x — 1 we find [TJ 22.4.2], 
(A.ll) 

(2 7 + n- l)(2 7 + n-2) 



(2 7 +l) 



T(2 7 + n) 
i!T(2 7 +l) 



for n > 2, with G^ 7) (l) = g[, 7) (1) = 1, where T(z) is the standard gamma function 
PP. Note that Gn \l) > for 7 > — 1/2 and that it decreases with increasing n if 
— l/2< 7 <l/2 but increases with n if 1/2 < 7 . 
Now (|A.5[) gives 



(A.12) 



d k G^ ] 



(x) 



2 fe r( 7 + fc + i) (7+fc) 



dx fe v ~ y r( 7 + i) 

which coupled with (|A.11[) gives 

2 fe - 1 r( 7 + fc)r(n + 2 7 + fc) 



(A.13) 



dx k 



(1) 



(n- fe)!r( 7 + l)r(2 7 + 2fc) ! 



n > 1. 



Gegenbauer polynomials correspond to Chebyshev polynomials of the 1st kind, T n (x), 
when 7 = 0, to Legendre P„(a;) for 7 = 1/2 and to Chebyshev of the 2nd kind, t/ n (a;), 
for 7 = 1. For the non standard normalization, 



(A 



14) G W(x) = ^ ) G^(x) = P n (x), <#>(*) = 



U n (x) 



Appendix B. Integrals and Asymptotics. As shown in section [3J the tau 
equations (13. 1| provide a matrix eigenproblem of the form fiAa — Ba. To reduce the 
coefficients A(k, I) and B(k, I) defined in (|3.2j) and (13. 3|) to the expressions Q3.4p and 
(|3~5]) . use (IA.4j) and (|A.5|) repeatedly 



(B.l) D 2 
1 

_ ~5 



2^2^(5/2) 



(X) 



-D 2 



(l + l)(l + 4)D \(l-x 2 )G\ 3 + /?(x) 



(l-x 2 ) 2 DG^(x) 
1 
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(i + l)(Z + 4)D (l-x^G<y 2 2) (x) 



^(/ + l)a + 2)(; + 3)a + 4)G^ 2 2) (x) 



(1/2) 



i <~ r Z+2 



Or). 



For the perturbation analysis described in section [3] we need the first order corrections 
to B(0, n - 4) and to B(l, n - 5). From equations ([3~T2]) and ([3~T3]) 



(B.2) 
(B.3) 



Bi(0,n-4) =lim 



Cn— 4 

e 

S 1 (l,n-5)=lim — 



P„_ 2 (a;) (l-x 2 ) 6 dx, 
x P„_ 3 (a;) (1 - x 2 ) c rfx. 
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To evaluate P n (x)(l~x 2 ) e dx to 0(e) for |e| < 1 and n even, use (IA.4I) for 7 = 1/2 
and integration by parts to derive 



(B.4) 



/ P n {x)(\-x 2 Ydx= 1 f D[(l-x 2 )DP n ](l-x 2 Ydx 
J -I n[n + 1 j J _i 

= -^j\p n (l-x 2 y(-2x)dx. 



This integral is if n is odd since P n (x) = (—l) n P n (—x). Since we have an e pre- 
factor we can now set e = in the integral, and do the remaining integral by parts to 
obtain 

"1 1, el 



/ P n (x)(l - x 2 Ydx — / xDP n dx 



-2e ' x 



( B -5) = — t / (£>(a:P„) - P„) 

n(n + 1) y_! 

^ :(P»(1) + Pn(-1)) 



n(n + 1) n(n + 1) 

for n even (0 for n odd as should be). 

For the integral in (|B.3|) , use the recurrence (|A.6|) for 7 = 1/2 to write (2n — 
5)xP n -3(x) — (n — 2)P„_ 2 (aO + {n— 3)P n ~i{x) and evaluate the resulting 2 integrals 
from the formula (IB. 51). Hence, for n even, 



AC 

(b.6) /;,>()./, i:- " : 



(B.7) Bi(l,n-5)~- 



(n-2)(n-l)' 



(n - 4) (71 - 1) ' 

Now for Ao(0, n — 4) and Ao(l, n — 5) and n even, we have 
(B.8) A {Q,n- A) =C n - 4 J D 2 P n _ 2 {x)dx = (n - 2)(n - 1)C„_ 4 , 

(B.9) A)(l,n-5) =C n -5 J xD 2 P n - 3 (x)dx = (n-4)(n~l)C n - 5 . 

Appendix C. Proof of theorem 15.51 

Consider 

00 

(C.l) /„(*; z)=J2 z k D k G^(x) + K (g^\x) - G%{x) 

k=0 

where z is a solution of f n (l; z) — and K is 
(C.2) K 



(0^,(1) -G%(1)) n + 2 



2(n + 7 )(Gl 7) (l)-G^ 2 (l)) "' 2(n + 7 + l)(n + 2 7 -l) 

where we have used (|A.11[) . Note that f n (l;z) — pn\z) defined in theorem 15.51 
Taking the x-derivative of f n (x; z) and using (|A. 10|1 . we find 



(C3) = ]T z k D k +^\x) 2K(n + 7 + l)G™i(*) 

k=0 
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Thus f n (x; z) satisfies the following differential equation 

(C.4) f n (x; z) - {l + K)G%\x) + KG%{x) = z^^- + z2K{n + 1+ l)G^ +1 {x)- 

Multiplying by (1 + x) < ^" t j^' 2 ' ) , integrating in the Gegenbauer norm and adding the 
complex conjugate we get 



1 d\f n \ 



(C.5) / ^p-(l+x)w(x)dx-(l + K) / {l+x)^-G^\x)w(x)dx + C.C. 



df* 



+ K (^J 1 (1 + x)^G%l 2 (x)w(x)dx + C.c)j =(z + z*) J 



i: 


dfn 




dx 



(1 + x)w(x)dx 



dfn n h) 



+ [z2K{n + 1 +l) / [l + x)^G™ x {x)w{x)dx + C.C 



where C.C. denotes the complex conjugate. To simplify (|C.5[) we need to compute 
four simple integrals 



1 d\f n 



_i dx 



-(1 + x)w(x)dx 



(C.6) 



J„ = / (l + x)^-G^(x)w(x)dx 



Jx 



{l + x)^G^> rl {x)w{x)dx 



i df n 

J 2 = J\l + x) d ^G%{x)w{x)dx 
Using integration by parts, the first integral becomes 



(C.7) h = j^ 



1 d\f n \ 



dx 



■(1 + x)w(x)dx 



1 \f n f(l-2 7 x)p^-dx. 
i 1 — x 



Therefore, the integral is negative for — 1/2 < 7 < 1/2 since for this range of param- 
eters 1 — 272; is positive. 

For the calculation of the three other integrals we are going to need the expression 

( C - 8 ) = S zkDk+lG n } (*) - 2K ^ + « + 

k=0 

(C.9) =V n - 2 (x; z) + 2( 7 + n - ^G^x) - 2K( 1 + n + ^G^x) 



(CIO) 



--V n -i(x; z) - 2K( 1 + n+ l)G^(x) 



where V n -2(x\ z) and V n -i(x; z) are polynomials of degree n— 2 and n— 1, respectively. 



GEGENBAUER TAU METHODS WITH AND WITHOUT SPURIOUS EIGENVALUES 19 



With the use of (|C.9|) and orthogonality of the Gegenbauer polynomials, we find 



(C.ll) J = J (l + x)^G%\x)w(x)dx 



2(71-1+7) / xG^l 1 {x)G^(x)w{x)dx-2K(n+l+j) / xG%L{x)G£\x)w{x)dx 



nj (G^\x)fw{x)dx-K{n + 1 + l)^- f {G^l 1 {x)fw{x)dx 

n T 7 J — 1 



^(7) _ g (n+7+l)(»+l) ^) 
n + 7 " +1 



where we have also used (|A.5[) and (|A.6[) . In the same way we compute 



(C.12) J 2 



i (l+^) T ^G^ 2 ( a ;) U ;(x)dx = 

(t) („\M-y) 



K / 2(7i + 7 + l)xG^' 1 {x)Gi;Ux)w{x)dx = -K(n + 2)h 



-1 



(7) 

n+2 



and 
(C.13) 

Ji = 



, d fn r Sl) 

dx 



(l+x)^G%l 1 (x)w(x)dx = -2K{n+i+l) I (l+x)(G { ^ 1 (x)) 2 w(x)dx 



- 2K(n + 7 + 1) / {G^Ux)fw{x)dx = -2K(n + 7 + lR+i 



(7) 



Thus, (|C.5[) transforms to 
(C.14) 

| 2 (1 - 27x)w;(a;) 



(1-z) 



- 2(1 + if) ( nh^ - /f (n + 7 + 1)(n+1) ^ 



n+l 



-2i<: 2 (n + 2)/7, 



(7) 

n+2 



f/.j; 



n + 7 

(1 + x)w(a;)da; - 4fcT 2 (n + 1 + 7) 2 ^ 1 



Our task is to show that the left hand side of the above expression is negative whereas 
the coefficient of the term z+z* on the right hand side is positive. A simple calculation 
shows that 

(C15) 77,^) - ^ + 7 + 1)^+1) ^ = nhi ) (» + !)(« + 2) { ) 
(77 + 7) " +1 " (n + 7)(77 + 2 7 - 1) 71+1 



7r 2- 1 - 2 TT(n + 27) 
-y 2 r 2 (7)77! (77 + 7) 



> 



7r2 -l-2 7r( ' n+27 ) 

7 2 r 2 (7)77,1(77 + 7) 



(n + 2)(n + 2 7 ) 
2(n-l + 2 7 )(7i+l + 7 ) 

(77 + 2)(77 + 1) 



(77 - 2) (277 + 1) 
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The last parenthesis is positive for n > 3. 

For the right hand side we use the notation in (|C.10[) to get 



(C.16) 



L 


dfn 




dx 



(1 + x)w(x)dx - 4K 2 (n + 1 + 7) 2 /ii+i 



(1 + x) (\V n -i{x; z)\ 2 - 2K(n + 1 + 7)^1 (x) z) + V* n _ x {x; z)) 

1 v 

+4K 2 {n + 1 + 7) 2 (Gi 7 i(a:)) 2 ) w(x)dx -4K 2 {n + l + lfh% = 
(1 + x)\V n -i(x; z)\ 2 w(x)dx + + + 4K 2 (n + 1 + 7) 2 ^+i - ^K 2 (n + 1 + ffh^ 



1 

(1 + x) \V n -i(x; z)\ 2 w{x)dx > 0. 



Thus (|C.14j) becomes 
(C.17) 



1 IAI 2 (1 - j ^ d X -2(l + K) (nW _ g (" + 7 + D(n + 1) e ) 
1 I 1 \ n-t-7 / 

- 2if 2 (n + 2)4+2 = (z + z*) / (l + a;)|P n _i(a;;2)| 2 «;(ar)da; 



and $t(z) < 0. 

For n — 2 we get the following characteristic polynomial 

(C.18) /| 7) (1; z) = H( 7 + i)( 3 z 2 + 3z + 1) 
whose zeros 

en m\ 1 V^i 1 , V3i 

(C19) z\ = . z 2 = 

\ ) 1 2 6 ' 2 6 

have negative real parts for any 7. 
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